Adaptive Grid Refinement for Discrete Tomography

Discrete tomography has proven itself as a powerful approach to image reconstruction from limited data. In recent years, algebraic reconstruction methods have been applied successfully to a range of experimental data sets. However, the computational cost of such reconstruction techniques currently prevents routine application to large data-sets. In this paper we investigate the use of adaptive refinement on QuadTree grids to reduce the number of pixels (or voxels) needed to represent an image. Such locally refined grids match well with the domain of discrete tomography as they are optimally suited for representing images containing large homogeneous regions. Reducing the number of pixels ultimately promises a reduct

Robust artefact reduction in tomography using Student’s t data fitting

Algebraic methods are popular for tomographic image reconstruction from limited data. These methods typically minimize the Euclidean norm of the residual of the corresponding linear equation system. The underlying assumption of this approach is that the noise has a Gaussian distribution. However, in cases where large outliers are present in the projection data, e.g., due to defective camera pixels, photon starvation from metal implants etc., the equation system is not consistent and the reconstruction will be fitted to these outliers, resulting in artefacts in the reconstruction. In this paper we use a penalty function for the residual that is based on the maximum likelihood estimate from the Student’s t distribution, which assigns a smaller penalty to outliers. No preprocessing is required to locate the outliers. We demonstrate the effectiveness of this approach on a 3D cone-beam simulated dataset for a series of perturbations in the projection data. Our results suggest that artefacts due to metal objects, de

Forty-fifth Annual Report for the Fiscal Year Ending June 30, 1968

Binary tomography is concerned with the recovery of binary images from a few of their projections (i.e., sums of the pixel values along various directions). To reconstruct an image from noisy projection data, one can pose it as a constrained least-squares problem. As the constraints are nonconvex, many approaches for solving it rely on either relaxing the constraints or heuristics. In this paper, we propose a novel convex formulation, based on the Lagrange dual of the constrained least-squares problem. The resulting problem is a generalized least absolute shrinkage and selection operator problem, which can be solved efficiently. It is a relaxation in the sense that it can only be guaranteed to give a feasible solution, not necessarily the optimal one. In exhaustive experiments on small images (2 × 2, 3 × 3, 4 × 4), we find, however, that if the problem has a unique solution, our dual approach finds it. In the case of multiple solutions, our approach finds the commonalities between the solutions. Further experiments on realistic numerical phantoms and an experiment on the X-ray dataset show that our method compares favorably to Total Variation and DART

Variable projection for non-smooth problems

Variable projection solves structured optimization problems by completely minimizing over a subset of the variables while iterating over the remaining variables. Over the last 30 years, the technique has been widely used, with empirical and theoretical results demonstrating both greater efficacy and greater stability compared to competing approaches. Classic examples have exploited closed-form projections and smoothness of the objective function. We extend the approach to problems that include non-smooth terms, and where the projection subproblems can only be solved inexactly by iterative methods. We propose an inexact adaptive algorithm for solving such problems and analyze its computational complexity. Finally, we show how the theory can be used to design methods for selected problems occurring frequently in machine-learning and inverse problems

Comparing RSVD and Krylov methods for linear inverse problems

In this work we address regularization parameter estimation for ill-posed linear inverse problems with an penalty. Regularization parameter selection is of utmost importance for all of inverse problems and estimating it generally relies on the experience of the practitioner. For regularization with an penalty there exist a lot of parameter selection methods that exploit the fact that the solution and the residual can be written in explicit form. Parameter selection methods are functionals that depend on the regularization parameter where the minimizer is the desired regularization parameter that should lead to a good solution. Evaluation of these parameter selection methods still requires solving the inverse problem multiple times. Efficient evaluation of the parameter selection methods can be done through model order reduction. Two popular model order reduction techniques are Lanczos based methods (a Krylov subspace method) and the Randomized Singular Value Decomposition (RSVD). In this work we compare the two approaches. We derive error bounds for the parameter selection methods using the RSVD. We compare the performance of the Lanczos process versus the performance of RSVD for efficient parameter selection. The RSVD algorithm we use i

A regularised total least squares approach for 1D inverse scattering

We study the inverse scattering problem for a Schrödinger operator related to a static wave operator with variable velocity, using the GLM (Gelfand–Levitan–Marchenko) integral equation. We assume to have noisy scattering data, and we derive a stability estimate for the error of the solution of the GLM integral equation by showing the invertibility of the GLM operator between suitable function spaces. To regularise the problem, we formulate a variational total least squares problem, and we show that, under certain regularity assumptions, the optimisation problem admits minimisers. Finally, we compute numerically the regularised solution of the GLM equation using the total least squares method in a discrete sense

Seismic wavefield redatuming with regularized multi-dimensional deconvolution

In seismic imaging the aim is to obtain an image of the subsurface using reflection data. The reflection data are generated using sound waves and the sources and receivers are placed at the surface. The target zone, for example an oil or gas reservoir, lies relatively deep in the subsurface below several layers. The area above the target zone is called the overburden. This overburden will have an imprint on the image. Wavefield redatuming is an approach that removes the imprint of the overburden on the image by creating so-called virtual sources and receivers above the target zone. The virtual sources are obtained by determining the impulse response, or Green's function, in the subsurface. The impulse response is obtained by deconvolving all up- and downgoing wavefields at the desired location. In this paper, we pose this deconvolution problem as a constrained least-squares problem. We describe the constraints that are involved in the deconvolution and show that they are associated with orthogonal projection operators. We show different optimization strategies to solve the constrained least-squares problem and provide an explicit relation between them, showing that they are in a sense equivalent. We show that the constrained least-squares problem remains ill-posed and that additional regularizati

A distributional Gelfand–Levitan–Marchenko equation for the Helmholtz scattering problem on the line

We study an inverse scattering problem for the Helmholtz equation on the whole line. The goal of this paper is to obtain a Gelfand–Levitan–Marchenko (GLM)-type equation for the Jost solution that corresponds to the 1D Helmholtz differential operator. We assume for simplicity that the refraction index is of compact support. Using the asymptotic behavior of the Jost solutions with respect to the wave-number, we derive a generalized Povzner–Levitan representation in the space of tempered distributions. Then, we apply the Fourier transform on the scattering relation that describes the solutions of the Helmholtz scattering problem and we derive a generalized GLM equation. Finally, we discuss the possible application of this new generalized GLM equation to the inverse medium problem

A data-driven approach to solving a 1D inverse scattering problem

In this paper, we extend a recently proposed approach for inverse scattering with Neumann boundary conditions [Druskin et al., Inverse Probl. 37, 075003 (2021)] to the 1D Schrödinger equation with impedance (Robin) boundary conditions. This method approaches inverse scattering in two steps: first, to extract a reduced order model (ROM) directly from the data and, subsequently, to extract the scattering potential from the ROM. We also propose a novel data-assimilation (DA) inversion method based on the ROM approach, thereby avoiding the need for a Lanczos-orthogonalization (LO) step. Furthermore, we present a detailed numerical study and A comparison of the accuracy and stability of the DA and LO methods

Full-waveform inversion with Mumford-Shah regularization

Full-waveform inversion (FWI) is a non-linear procedure to estimate subsurface rock parameters from surface measurements of induced seismic waves. This procedure is ill-posed in nature and hence, requires regularization to enhance some structure depending on the prior information. Recently, Total-Variation (TV) regularization has gained popularity due to its ability to produce blocky structures. Contrary to this, the earth behaves more like a piecewise smooth function. TV regularization fails to enforce this prior information into FWI. We propose a Mumford-Shah functional to incorporate the piecewise smooth spatial structure in the FWI procedure. The resulting optimization problem is solved by a splitting method. We show the improvement in results against TV regularization on two synthetic camembert examples